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Abstract 

The transcription factor p53 is an important regulator of cell fate. Mutations in 
p53 gene are associated with many cancers. In response to signals such as DNA dam- 
age, p53 controls the transcription of a series of genes that cause cell cycle arrest 
during which DNA damage is repaired, or triggers programmed cell death that elimi- 
nates possibly cancerous cells wherein DNA damage might have remained unrepaired. 
Previous experiments showed oscillations in p53 level in response to DNA damage, but 
the mechanism of oscillation remained unclear. Here we examine a model where the 
concentrations of p53 isoforms are regulated by Mdm22, Arf, Siah, and /3-catenin. The 
extent of DNA damage is signalled through the switch-like activity of a DNA damage 
sensor, the DNA-dependent protein kinase Atm. This switch is responsible for ini- 
tiating and terminating oscillations in p53 concentration. The strength of the DNA 
damage signal modulates the number of oscillatory waves of p53 and Mdm22 but not 
the frequency or amplitude of oscillations-a result that recapitulates experimental find- 
ings. A critical finding was that the phosphorylated form of Nbsll, a member of the 

DNA break repair complex Mrell-RadSO-Nbsll (MRN), must augment the activity of 
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Atm kinase. While there is in vitro support for this assumption, this activity appears 
essential for p53 dynamics. The model provides several predictions concerning, among 
others, degradation of the phosphorylated form of p53, the rate of DNA repair as a 
function of DNA damage, the sensitivity of p53 oscillation to transcription rates of 
SIAH, /3-CATENIN and ARF, and the hysteretic behavior of active Atm kinase levels 
with respect to the DNA damage signal. 

Keywords: DNA damage, p53, Mdm2, oscillation, nonlinear dynamics, cancer, systems 
ology, computational model 



The p53 protein is a transcription factor that is present in most higher eukaryotes. It 
regulates the expression of a series of genes that mediates cell cycle arrest and apoptosis in 
response to DNA damage and extracellular signals p53 gene is often mutated in tumor 
cells 12 0] ■ It is necessary for the elimination of potentially cancerous cells by triggering cell 
death |3j and is thought to be a major tumor suppressor gene p53 protein is subject 

to a series of post translational modifications that modulate its kinase activity jl] . A simple 
view is that p53 is present in at least two major forms, an unphosphorylated form that is 
thought to be inactive as a transcription factor, which is rapidly turned over by proteolysis 
and a more stable phosphorylated form. Phosphorylation occurs at a number of amino acid 
residues on the protein by the activity of protein kinases, and the phosphate groups are in 
turn removed by the activity of several phosphatases that are dynamically associated with 
p53. 

A phosphorylated form of p53, p53-P, is the active transcription factor that binds to a 
number of target DNA sites that regulate the expression of several genes involved in cell cycle 
checkpoint arrest and programmed cell death, respectively [7]. Additionally, p53-P induces 
transcription of the gene p53AIPl, whose product in the cytoplasm leads to the release of 
cytochrome c protein from mitochondria, which in turn activates programmed cell death 
directly by activating the cytoplasmic protease cascade [H]. 

Ionizing radiation causes DNA double strand breaks that triggers the binding of the sen- 
sor protein kinase Atm [9 to the nascent DNA ends lOj. DNA binding triggers the kinase 
activity of Atm, which then phosphorylates preexisting p53 proteins to p53-P E] • Un- 
phosphorylated p53 is the substrate of Mdm2, a ubiquitin E3 ligase [12]. Once ubiquitinated, 
p53-P protein is proteolytically degraded by the proteosome complex. By contrast, p53-P is 
more resistant to ubiquitination, therefore is more stable. Hence Atm activation leads to the 
accumulation of p53-P. Nuclear localized p53-P binds to the regulatory sequences of several 
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genes among which MDM2 is one. This binding increases the transcription of MDM2 gene. 
Thus an increase in p53-P concentration in the nucleus increases the synthesis of Mdm2 RNA 
and protein. The increase in Mdm2 protein concentration is expected to destabihze any 
remaining unphosphorylated p53 by ubiquitination and proteasome-mediated degradation. 
Using this negative feedback, p53 regulates its own concentration in the cell |3]. 

DNA double strand breaks induced by ionizing radiation are thought to cause certain topo- 
logical changes in the chromatin structure such that the inactive Atm kinase homodimer 
molecules associated with the chromatin become activated ^2)- One of the two Atm kinase 
subunits phosphorylates the other at the serine-1081 residue leading to reciprocal phos- 
phorylation of the two Atm kinase subunits. The dimer is thought to dissociate and the 
phosphorylated monomers have three central functions: First, phosphorylation of nuclear- 
localized p53 to p53-P; second, phosphorylation of Brcal and Brca2 proteins to activate 
cell-cycle checkpoint arrest and blockage of DNA replication fork movement; and third, re- 
cruitment of the MRN (Mrell-Rad50-Nbsl) DNA repair complex to the nascent DNA ends 
by protein-protein interaction while simultaneously activating the function of the MRN com- 
plex by phosphorylating Nbsl to Nbsl-P ^1^]. The activated MRN complex begins to 
repair the broken DNA either by nonhomologous end joining or by recombinational repair 
using the unbroken homolog as a template. p53-P directly recruits RAD51 protein at the 
MRN complex by protein-protein interactions 16J where RAD51 constrains the recombina- 
tional repair of DNA by the non-crossover mode. This mode of DNA repair protects the 
genome from somatic aneuploidy that is often an early event during cancer formation [17]. 
Yet another effect of p53-P is that it slowly accumulates in the cytoplasm where it causes 
the mitochondrial membrane to leak cytochrome c protein fragments [8 . Accumulating 
cytochrome c in the cytoplasm is ultimately responsible for inactivating the anti-apoptotic 
proteins Bcl2/Bax, which leads to apoptosis. However, the relative rates of cytoplasmic 
p53-P accumulation and its nuclear effects in DNA repair determine the fate of the cell: 
repair of DNA damage or death. Thus, the levels of p53 protein and its phosphorylated form 
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are crucial determinants of cell fate. 

Recently Lahav et al. ^] measured intracellular concentration of total p53 and Mdm2 pro- 
teins by fusing their coding sequences to fluorescent reporter domains. Examination of single 
mouse cells following treatment with ionizing radiation that produces DNA double stranded 
breaks (DSB) revealed that p53 and Mdm2 protein concentrations in singles cells oscillate 
in response to DNA damage. The period and the amplitude of these oscillations remained 
constant over a range of radiation doses, but the number of pulses was proportional to the 
radiation dose. Lahav et al. ^Hl proposed that the system behaves as a "digital oscillator". 

In response, Tyson [19J proposed that following DNA damage the steady state p53 concen- 
tration passes through a Hopf bifurcation fl^ m\ and begins to oscillate. Once the damage 
is corrected the system is pulled back from the oscillatory region to its original steady state. 
This idea was further expanded in an explicit model that assumed regulation of p53 degra- 
dation by Mdm2 and p53 as a transcription factor for Mdm2 RNA synthesis [22]. A crucial 
aspect of this model is a positive feedback loop by which p53 inhibits the nuclear localization 
of Mdm2 into the nucleus [221 , thereby increasing the p53 levels. 

Another approach to modeling the p53 dynamics makes explicit use of delays in the sys- 
tem corresponding to the time it takes for transcription and translation of the Mdm2 pro- 
tein 23 ^- Oscillations are possible by including a sufficiently long delay, as for example in 
the case of a hypothetical intermediate step introduced by Alon et al. j2S]. 

While previous modeling approaches explored the conditions that could in principle gener- 
ate oscillations in p53 and Mdm2 levels, these models did not specifically address explicit 
molecular mechanisms of DNA break repair or control of p53 activity. Here, in a significant 
departure from previous studies, we explicitly model changes in unphosphorylated and phos- 
phorylated versions of p53 in response to DNA damage. We ignore the nuclear transport 
of p53 because the time-scale of nuclear localization (minutes) is fast compared to that of 

the observed oscillatory dynamics (hours). We include in the model a specific pathway, that 
iTiana, G., Jensen, M.H., & Sneppen, K., (2002) arXw:cond-mat/0207236i )l 
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includes the regulation of Mdm2 through p53-P using Arf(pl9), /3-catenin, and Siah 3]. We 
also explicitly model the activation of the Atm kinase, which phosphorylates p53 as a switch. 

The model reproduces the experimentally observed oscillatory dynamics in p53 and Mdm2 pro- 
tein concentrations and generates a range of novel molecular predictions that can be tested by 
future experiments. These predictions include: the expected degradation profile of p53-P, the 
rate at which DNA is repaired is proportional to the degree of DNA damage, the sensitivity 
of p53-P oscillations to the transcription rates of several genes including SIAH, /?-CATENIN 
and ARF, and, finally, specific predictions on the hysteretic behavior of the Atm kinase 
switch. 

Methods 

All simulations were carried out using the Systems Biology Workbench (SBW/BioSPlCE) 

tools [21]: the network designer, JDesigner, the simulation engine Jarnac ^7\. Bifurcation 

diagrams were computed using SBW with an interface to MATLAB j2H], and a bifurcation 

discovery tool [21] ■ Bifurcation plots were also computed and cross checked using OscillS ^, 

an interactive bifurcation software package which is linked to AUTO jSOj, and SBW [26j. 

In all our simulations the species concentrations are regarded as dimensionless, whereas 

the kinetic constants have dimensions of inverse time, with dimensionless Michaelis-Menten 

constants. All models are available as standard SBML or Jarnac scripts. In the sensitivity 

analysis of a network, the control coefficients Cj are defined as, Cj = where Si, pj are 

the species, and the parameters respectively, 
^http: / / sourceforge.net /projects /oscillS 
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Results 



The Model 

We modularize the network j^I] into three distinct dynamically different pieces: a switch 
describing the Atm kinase, a subnetwork that describes the DNA damage repair, and an 
oscillator that includes p53 and Mdm2 (Fig 1). DNA damage turns on the Atm kinase 
switch which then turns on the oscillator. The oscillations activate and maintain the repair 
process. As the damage is corrected, the switch is shut off; this in turn switches off the 
oscillator (see Fig 6). 

Module A: DNA Repair 

We assume that when DNA is damaged a signal Sn, which is proportional to the damage, is 
generated. Repair of DNA breaks decreases Sd- In reality, Sd is equivalent to the number of 
nascent double stranded DNA ends created due to the interaction of ionizing radiation with 
DNA. In response to the DNA damage, Atm is phosphorylated to Atm-P, and the activated 
Atm-P kinase in turn phosphorylates a host of target proteins, including Nbsl, Brcal, Brca2 
and p53. Brcal and Brca2 block DNA replication and cell cycle. The phosphorylated form 
Nbsl-P recruits Mrell and RadSO, and the complex of three proteins (the MRN complex) 
binds the broken DNA to initiate repair-recombination. Hence the repair process begins 
once concentrations of Atm-P kinase, Nbsl kinase and p53-P begin to grow. The rate of 
repair is thus described by 

^ = -a^[p53-P] [Nbsl-P] [Atm-P] (1) 

where ac = 2 x 10^'^, is a rate constant of the decay in damage signal, which can be calculated 
by measuring the change in the number of DNA ends per genome. The bracketed variables 
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refer to the concentration of the various species, in dimensionless form. For simphcity we 
assume that all double strand break (DSB) repair is dependent on p53; this is clearly untrue 
because more error prone p53 independent repair pathways exist. 

Module B: Atm Switch 

Atm proteins are normally sequestered in the form of dimers/multimers, and are thought to 
undergo autophosphorylation jT2lE21- A class of phosphatases, PP-2A, are known to de- 
phosphorylate Atm-P [231 • When a DSB occurs, two events happen very quickly. The first is 
that the phosphatases dissociate from Atm multimers and this results in an increased phos- 
phorylation rate, and the second is that the Atm dimers dissociate to activated monomers. 
It has been hypothesized that "turning off" the suppressive effect of the phosphatases on 
the autophosphorylation of Atm allows rapid signal transduction of the DNA damage sig- 
nal [Sni- Upon DNA damage, Atm-P accumulates rapidly in the nucleus jT2j. We assume 
that the effect of ionizing radiation is to reduce the ability of the phosphatases to bind to 
Atm. This reduction in the phosphatase activity following DNA damage could be thought 
to occur due to competition among Atm-P and the signal Sd to bind to the phosphatases. 
However, we do not explicitly model the competitive inhibition but use a mathematical rate 
law that is derived from such an assumption. It is known that Atm-P phosphorylates Nbsl, 
which is part of the MRN complex, and that the MRN complex further increases the activ- 
ity of Atm ^j. We take this as a cue to modelling the interaction between Atm and 
Nbsl as a positive feedback system. As such, active Atm (i.e. phosphorylated Atm) will 
activate Nbsl in the MRN complex (through phosphorylation), leading to further activity 
in Atm. This explains why Atm-P levels increase rapidly in the cell upon DNA damage. 
Thus, Nbsl-P positively regulates (thick dotted line in Fig. 1, Module B) activity of Atm by 
causing its phosphorylation. Alternatively, one can regard this effect to be due to a further 
activation of Atm-P by Nbsl-P; for simplicity we do not introduce a new variable for the 
"active" form of Atm-P. We do not know the exact process by which this may occur, but 
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find that such a process is required by the model to be able to predict the experimental 
observations of the p53 dynamics (see later). 

Referring to Fig. 1, Module B, the equations (see supplementary materials for reaction 
scheme) for Atm, and Nbsl are, 



fi[Atm-P] 
di 

ci[Nbsl-P] 
Jt 



[Atm][Nbsl-P] 
ki + [Atm] 

[Nbsl] [Atm-P] 
h + [Nbsl] 



+ a 



+ f3 



[Atm-P] 
k2 + [Atm-P] 
[Nbsl-P] 
+ [Nbsl-P] 



D 



(2) 



where the assumption of the positive feedback (marked as the dotted line in Fig. 1 Module B) 
is given by the first term in the rate law for Atm-P. a, P, are basal rates of phosphorylation 
due to stochastic fluctuations in Atm and Nbsl kinase activities, respectively. In the above 
equations, the rate constants are in units of inverse time, is a scaling parameter, and 
the Michaelis-Menten constants [ki, k2, k^, k4^, ko) are in dimensionless units. The parameter 
values for these equations are reported in supplementary materials. 

Using the above equations the steady state values of Atm-P reached by the system as a 
function of the DNA damage, Sn, was computed in Fig. 2 A. The system shows bistability, 
i.e. two possible steady states for a range of Sd- For large values of damage {Sd > 2.7), 
Atm-P = 1, and the switch is in the ON state. As the damage is corrected, and drops below 
~ 0.08, the switch is turned OFF. The plot for Nbsl-P is very similar, and this is due to the 
nature of the positive feedbacks that Atm and Nbsl have on each other. There are several 
examples where positive feedback can give rise to bistable, switch-like, behavior |Mll35| l3H]. 
The system can therefore move between two states, a low state defined by low values of 
Atm and Nbsl, and a high state defined by high values of Atm-P and Nbsl-P. As shown in 
Fig. 2B, as the damage decreases below ^ 0.08, the switch is turned to the OFF state, and 
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Atm-P levels drop off rapidly. 

In the present model the positive feedback between Atm-P and Nbsl-P has significant effects 
on the dynamics of the integrated system. Neglecting the positive feedback still gives a quick 
response in Atm-P levels to DNA damage, but as the DNA is repaired, the Atm-P levels do 
not fall off abruptly, but gradually decrease (supplementary materials). For the integrated 
model to produce the observed jJHI results that describe the experimental observations it 
is necessary to have a switch-like behavior in Atm-P levels that in turn requires the stated 
positive feedback. 

Module C: p53 Oscillator: 

To determine whether the known regulatory network involving p53 and Mdm2 is sufficient 
to reproduce their observed dynamics, we proceed as follows. The central oscillator is com- 
posed of the negative feedback of Mdm2 on p53 stability through ubiquitination, and the 
dependence of Mdm2 expression on p53-P. The latter is produced by Atm-P in response 
to the DNA damage. The dephosphorylation of p53-P is described by a Michaelis-Menten 
scheme, and is assumed to occur at a very small rate compared to the phosphorylation rate. 
We assume small first order degradation rates for both p53 and p53-P. p53-P is a more 
stable form of p53, and hence its rate of degradation is smaller than the degradation rate 
for p53. We model the degradation of p53-P as a function of Mdm2; this is done to match 
qualitatively the simulation results of our model to the experimental results (see below). 
Furthermore, p53 indirectly regulates the amount of Mdm2 that is available for binding to 
p53 through the former's transcription factor activity on the SIAH gene jl]. The Siah pro- 
tein is a ubiquitin ligase, and enhances the degradation of /5-catenin. /3-catenin, which is 
produced constitutively, regulates the transcription of the Arf gene. The Arf protein binds 
to Mdm2, sequestering the latter from its substrates and also promoting the latter 's degrada- 
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tion such that increasing concentration of Arf should increase the p53 level ^. The following 
differential equations describe the dynamics of module C in Fig. 1, and were generated 
automatically by SBW from SBML fl^ . 



rf[p53] 
dt 



d[p53-P] 
Jt 



d[Mdm2-p53] 
di 

(i[Mdm2-p53-P] 
di 



[Atm-Pl [p53l 
Co - —, — , r-^, - con + Ph 



[p53-P] 



kp + [p53] 
A;ip[p53][Mdm2] - A;2p[p53] 
[Atm-P][p53] 



Ph 



kp + [p53] 

[p53-P] 



{kdep + [p53-F]) 

+ con - A;p3[Mdm2][p53-P] 



(3) 



A:d[p53-P], 



{kdep + [p53-F]) 
A;ip[p53][Mdm2] - A;ci[Mdm2-p53], 

A;p3[p53-P][Mdm2] - A;c2 [Mdm2-p53-P] , 



rf[Mdm2] 
di 



fi[Mdm2-Arf] 

Jt 
d[Arf] 

dt 
d[Siah] 

dr~ 

(i[/3-catenin 
Jt 

(i[/3-catenin-Siah] 
Jt 



Pi , ^^^f"rjpi4 - ^ip[Mdm2][p53] - ^4,, [Mdm2] [Arf] - A;5[Mdm2] 
km + [p53-PJ'' 

A;p3[Mdm2][p53-P] + A;c2[Mdm2-p53-P] + A;ci [Mdm2-p53] , 
k^a [Mdm2] [Arf] - k^ [Mdm2-Arfl , 

Po[/5-catenin] - fc4a[Mdm2][Arf] - A;7[Arf] + fce [Mdm2-Arf] , 
[p53-P]^ 



P3 



+ A;c3[/?-catenin-Siah] — A;c4[/9-catenin] [Siah] — ^^^[Siah], 



A;, + [p53-P]4 
Pip — A;8[/5-catenin] — A;c4[/9-catenin][Siah], 

A;c4[/3-catenin][Siah] — A;c3[/3-catenin-Siah] 



The parameter values for these equations are given in supplementary materials. The tran- 
scriptionally active p53-P proteins form a tetramer jHZj, and hence the transcription rates 

^Sionov , R.V., Hayon, I.L., fc Haupt,Y., (2000 -2005) T he regulatio n of p53 growth 
supression |http://www.ncbi.nlin.nih. gov/books/bv.fcgi?rid=eurekah. chapter. 11736| 
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for MDM2 and ARF genes are described by a Hill coefficient of 4. The concentration of 
Atm-P, which determines the amount of phosphorylated p53, is treated as a parameter in 
these equations. 

In Fig. 3A, we plot the steady state concentrations of p53-P as a function of the Atm-P con- 
centration. p53-P exhibits oscillations for Atm-P > 0.22, i.e. in a region beyond the 
Hopf bifurcation. The kinetics for p53-P and Mdm2 concentrations, plotted in Fig. 3B for 
Atm-P = 1, show that they are out of phase with each other. p53-P is assumed to decay 
by two mechanisms: a ffist order decay that is independent of Mdm2 concentration, and 
another that is dependent on Mdm2 concentration. The second decay rate was assumed 
because without this critical assumption, the amplitudes of p53-P oscillations do not fall to 
low levels (supplementary material) as was observed in the experiment jl8j . 

A sensitivity analysis shows that if the transcription rate of MDM2 is decreased, then this 
can be balanced by a slight increase in the transcription rate of SIAH, and vice versa. The 
decrease in SIAH transcription leads to a release of Mdm2 that was previously sequestered by 
Arf, and this counteracts the decrease in Mdm2. The control coefficients, C^^'^ ^ , C^f'^ ^ 
are both negative, and are approximately of the same order of magnitude (see Supplementary 
materials). Increasing either of these two transcription rates leads to a decrease in p53-P, 
and hence they are mutually compensatory. 



Integrated Model 

The three modules are now combined into one integrated model. Eqs. 1, 2 and 3 are 
simulated with appropriate initial conditions to obtain the concentrations of the various 
species as a function of time. The initial conditions are the steady state concentrations of 
p53, p53-P, Mdm2, Arf, /3-catenin, Siah, and the complexes. These steady state values are 
obtained by solving the rate equations for So = 0. The DNA breaks are assumed to occur at 
t = 0. Although the DNA breaks occur over a finite time interval, we assume that this is short 



12 



compared to the period over which repair occurs. The DNA damage turns the Atm switch 
ON, which rapidly makes Atm-P = 1, and Nbsl-P = 1 (see Fig. 2B). Since Atm-P > 0.22 
(from Fig. 3A), the system crosses the Hopf Bifurcation point, and enters the oscillatory 
regime. p53-P begins to oscillate, and p53-P, Atm-P and Nbsl-P (through their respective 
actions not explicitly modelled here) begin to correct the damage (Eq. 1). As DNA repair 
proceeds. So decreases. At a certain threshold, Atm-P falls sharply back to its original value 
of zero (similarly for Nbsl-P). The Atm switch is shut OFF, and the system exits out of the 
oscillatory zone back to its original stable state. In Fig. 4A, we display the concentration 
of Atm-P, which rapidly rises to ~ 1, remains at this value for a while, and finally rapidly 
declines as So decreases below a certain threshold. Also shown is the dynamics of p53 which 
falls from its steady state value close to ~ 2, and oscillates with a small amplitude, but 
eventually rises to its original steady state value, as the repair is completed. All the other 
species described in Fig. 1 Module C also oscillate (i.e. Siah, /3-catenin, Arf etc). In Figs. 4 
(B-D), we plot the time series values of S^, p53-P and Mdm2 for different initial values of 
DNA damage. 

Note that the number of pulses increases with damage; however, the shape and frequency of 
the pulses do not change. This is consistent with the experimental observations that describe 
the "digital" behavior of p53 oscillations JH] • On DNA damage the system is moved into the 
oscillatory zone rapidly, by the Atm-P switch. Once inside the oscillatory zone, the switch 
is ON at a fixed value of Atm = 1, which fixes the amplitude of the oscillations. Similarly 
when the switch turns OFF the system moves out of the zone rapidly. Hence even though 
in this case the system exhibits a supercritical Hopf bifurcation, for which the amplitude 
of oscillations vary slowly with the value of Atm-P, the rapid movements into and out of 
the oscillatory regions give rise to fixed amplitude and frequency. An alternative model, as 
proposed in [22j is to employ a subcritical Hopf bifurcation which has the natural property 
to generate large amplitude oscillations with only small parameter changes. 
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Discussion 



The model presented here describes the dynamics of a gene-regulatory network involved in 
DNA damage repair. We describe the model as three interacting modules. The p53, Mdm2, 
Siah, /3-catenin, and Arf together represent the oscillator. The Atm, Nbsl dynamics presents 
a switch. The repair process is active as long as p53-P, Atm-P and Nbsl-P are present. By 
modularizing the network we identified the key dynamical assertions (Fig. ID). The Atm, 
Nbsl switch is activated by DNA damage; the switch turns on Atm-P and p53 is phos- 
phorylatcd; p53-P then begins to oscillate due to its interaction with Mdm2 (also through 
Arf, Siah, and /3-catenin). The damage begins to be repaired because p53-P, Atm-P and 
Nbsl-P are abundant. The number of pulses of p53-P is proportional to the damage but the 
frequency and amplitude of the pulses remain the same. As Sa decreases the Atm-P activ- 
ity is switched OFF and the system comes back to its steady state. Our model therefore 
supports a "digital" oscillatory behavior. 

An important point is the sensitivity of the dynamics of the integrated model to parameter 
values. In Fig. 5 we plot on a logarithmic scale the range of parameter values for which the 
model exhibits qualitatively similar behavior. 

There are a few critical parameters that determine the period and the amplitude of p53-P os- 
cillations. Specifically, the rates of transcription of SIAH and ARF have significant effects 
on the period. This is intuitively obvious because these rates determine the rate at which 
Mdm2 is released from sequestration due to binding with the Arf protein. p53-P amplitude 
can be increased by changing the Michaelis-Menten constant that determines the rate of 
phosphorylation (i.e kp). Similarly it can be decreased by increasing the transcription rate 
of MDM2 gene. These observations suggest simple tests. For example, partial functional 
knockout by siRNA could reduce the mRNA from any of the genes that form part of the 
p53 oscillator circuit. The model predicts that p53 oscillations would be affected by these 
partial knockouts. 
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The degradation of p53-P is a largely unknown process. If we only assume a simple first 
order rate of degradation proportional to p53-P the results of the simulations are inconsistent 
with experimental observation [TS]. In particular the troughs in the oscillations of p53-P do 
not reflect the low levels observed in the experiment. To accommodate this discrepancy, 
we proposed an additional degradation route which is a function of both p53-P and Mdm2. 
Although the mechanism assumed here is similar to ubiquitination, there may be another 
mechanism though an intermediate gene product. With this addition to the model, the simu- 
lations faithfully reproduce the experimental observations. Experiments to test this alterna- 
tive degradation routes should constitute removing Mdm2 and measuring p53-P degradation 
using antibodies specific to this isoform. 

A crucial aspect of the model is the positive feedback regulation between Atm-P and Nbsl-P. 
This leads to a bistable switch-like model for the Atm activity. Here we have not explic- 
itly modelled the mechanism of Atm activation by DSB, which may or may not require 
MRN. However, the positive feedback between Atm and Nbsl is essential. Most recent in 
vivo experiments provide additional support for the positive feedback activation of Atm by 
Nbsl [HH]- To verify the Atm switch behavior, it is necessary to determine the two damage 
thresholds that occur in the bifurcation plot (Fig. 2A): The first is the minimum damage 
required to throw the Atm switch ON {So — 2.7), and the second is the damage level at 
which the Atm switch turns OFF {Sd — 0.08). The first threshold can be measured if cells 
are subjected to increasing doses of irradiation, and the damage level that turns Atm-P level 
ON is determined. Once the Atm switch is ON, we wait for the switch to turn OFF and 
at that point measure the existing DNA damage. The model predicts that the time taken 
to correct a larger amount of damage will be almost as long as that for a small amount 
of damage. This is because for high Sd, Atm-P levels are high and the repair progresses 
relatively rapidly, whereas for low Sd, the Atm switch is OFF and so correction must occur 
through basal amounts of Atm-P. 

Transgenic mice with extra copies of p53 genes are resistant to tumor formation and age 
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normally, suggesting that increased doses of p53 gene may allow normal regulation pPj. By 
contrast, mice with a dominant negative form of p53 gene, whose product is always ON but 
is not subject to Mdm2 action, are radiation resistant but age faster. We predict that in 
these latter mice, p53 protein levels do not oscillate [^U] . 

We note (Fig. 6, curve B) that by adjusting a and (i (i.e., increasing the basal phosphorylation 
rates) it is possible to generate high levels of Atm-P even for a low Sd- By this analysis, 
if the basal rates of phosphorylation of Atm could be increased then the time taken to 
correct a small amount of damage could be decreased. This is testable by decreasing the 
concentrations of phosphatases through the use of siRNA. Furthermore, if is decreased, 
which corresponds to inefficient phosphatase binding to Atm, the steady state values of 
Atm-P for small Sd remain at a high level (Fig. 6, curve A). Therefore it would be 
instructive to examine the possible existence of cell lines where the levels of p53 naturally 
oscillate even in the absence of much DNA damage. We predict that, if such cell lines exist, 
mutations in these lines should affect the binding of phosphatases to Atm. 

In this work we have neglected the effects of stochastic fluctuations in chemical species by 
assuming that the concentrations of all species are large. In a future publication we will 
examine the effects of stochastic noise on small concentration ranges. 
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Supplementary Materials 

In this supplement we describe reaction schemes, parameter values, and supplementary fig- 
ures referred to in the main manuscript. 

Atm Switch Dynamics 

The reaction scheme, from which Eq 2, in the main text, are derived is described by, 

Atm Atm-P 
Atm-P Atm 
Nbsl Nbsl-P 
Nbsl-P ^ Nbsl 

Reaction Scheme for the Atm Switch 

The parameter values for these equations are reported in Table. 1. Each of Atm and Nbsl are 
conserved cycles, i.e. Atm-|-Atm-P = 1, and Nbsl-|-Nbsl-P = 1 



h 


k2 


^3 




ko 


a 




as 


0.01 


1 


0.1 


0.01 


0.1 


lE-4 


lE-4 


0.3 



Table 1: Parameters used for the Atm switch. 



In the model for the Atm switch, we assume a positive feedback, by allowing Nbsl-P to pos- 
itively regulate the phosphorylation of Atm to Atm-P. Here we show though simulation, the 
effect of neglecting the positive feedback, on the dynamics of Atm-P. Although Atm-P levels 
rise rapidly with DNA damage, the decline with subsequent DNA repair, is slow, and hence 



[Atm] [Nbsl-P] 
fci+[Atm] ^ 

[Atm-P] / a 
fe2+[Atm-P] Vfco+[6'i5] 

[Nbsl] [Atm-P] , o 
fc4+[Nbsl] 

[Nbsl-P] 
fc3+[N bsl-P] 
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this would not explain the switch-like behavior of Atm-P. As is explained in the main text, 
the rapid drop in Atm-P levels is required to turn OFF the p53 oscillator. 

If in the above reaction scheme, we leave out the positive feedback, in the activation of 

Atm-P by using Atm — > Atm-P , with the rate law, , , + «, then a simulation of the 

fcip+[/vtmj 

switch model leads to the results shown in Fig. 7. The dynamics of Atm-P does not rapidly 
decline as the damage decreases. The latter effect is required to turn off the p53 oscillator. 

In the main text in Fig. 6, we show plots of the steady state values of Atm-P, as a function 
of Sd, for two cases, case A: Modified phosphatase activity, with ko — 1. case B: Increase 
in the basal phosphorylation rates, a = 0.01,(3 = 0.01. The rest of parameters remain the 
same for both cases (as given in Table 1). 

p53 Oscillator Dynamics 

The equations for the p53 oscillator model given in Eq. 3, in the main text are obtained 
from the following reaction scheme. 
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source — > p53 

p53 ^ p53-P 

p53-P ^ p53 

p53-P sink 

source Mdm2 

Mdm2 + Arf ^Arf-Mdm2 

Arf-Mdm2^ sink + Arf 

p53 + Mdm2 ^p53-Mdm2 

p53-Mdm2^ Mdm2 + sink 

p53 — > sink 

source — > /3-catenin 

source — > Siah 

Siah + /3-catenin /3-catenin-Siah 

/?-catenin-Siah Siah + sink 

Arf sink 

source — Arf 

Siah sink 

/3-catenin — > sink 

p53-P + Mdm2 ^ p53-PMdm2 

p53-PMdm2 ^ Mdm2 + sink 

Mdm2 — > sinA; 



Co 

p [p53-P] 

(fed.j,+[p53-P]) 

fcrf[p53-P] 

pi [p53-P]^ 
(fc„+[p53-P]4) 

/c4a[Mdm2][Ar^ 

/c6[Arf-Mdm2] 

A;ip[p53][Mdm2] 

fcci[p53-Mdm2] 

fc2p[p53] 



[p53-P]^ 



(fc,+[p53-P]4) 
A;c4[Siah] [/3-catenin] 

/Cc3 [/3-catenin-Siah] 

fc7[Ar^ 

P6[/3-catenin] 

A;si[Siah] 

A;8[/3-catenin] 

fcp3[p53-P][Mdm2] 

A;c2[p53-P][Mdm2] 

A;5[Mdm2] 



where the parameter values are given in Table 2, , 



kp3 


kp 


con 


Ph 


kdep 


kd 


kip 


k2p 






Pi 




0.005 


0.4 


lE-4 


lE-6 


1 


0.006 


1 


0.5 


1 


1 


0.4 


1 


h 


kia 


kg 


ke 


Po 








kfA 


Pi 




kgi 


lE-4 


15 


lE+5 


2 


0.01 


0.3 


0.1 


2 


3 


1 


0.2 


0.02 



Table 2: Parameters used in the equations. 



and where the constitutive rate of production of p53 is Cq = 0.8. In Fig. 8A, we plot 
the bifurcation plot for the steady states of p53-P, as a function of Atm-P. There are 
two supercritical Hopf bifurcations, between which the system exhibits oscillations. Note 
however, that the minimum value of the p53-P oscillations do not tend to a low value, as is 
experimentally observed. The time series of the concentrations for p53, Mdm2 are displayed 
in Fig. SB, for Atm-P = 1.0. 

Sensitivity Analysis for module C 

In a sensitivity analysis of a network, the control coefficients Cj are typically computed, 

where Si, pj are the species, and the parameters respectively. Hence Cj expresses a fractional 
change in a species value, as one of the parameters is changed. The coefficient could be 
positive, or negative, depending on whether a small increase in a parameter tends to increase 
or decrease a species level. The sensitivity is obtained by linearizing the system about a 
steady state. In the case of the oscillator module, the system oscillates about an unstable 
point. We fix the value of Atm-P = 1, for which the steady state can be computed. The 
computation is done using the sensitivity features in JDesigner. In the table below, we 
compute the control coefficients for a few important parameters for the oscillator network 
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(we display only those control coefficients that are either reasonably large, or are important 
for regulation) , module C, which describe the regulation that occurs. . 



^p53-P 


^p53-P 
^Atm-P 


^p53-P 


^p53-P 


^p53-P 


0.16 


0.15 


-0.26 


-0.24 


0.26 



Table 3: Control coefficients for some of the important regulatory steps in module C. 



For example we see that Ceo is positive, which is intuitively obvious, since the parameter 
Co determines the production of p53. The control coefficients Cp^'^ ^, Cp^'^ ^ , are approx- 
imately of same order, and of the same sign. An increase in Pi,Ps represents an increase 
in the transcription of MDM2, SIAH genes, both of which lead to a reduction in p53-P. 
Hence a reduction in either transcription rate, can be compensated by an increase in the 
other. C^^"^ ^ is positive due to the following: an increase in Pb, leads to an increase in Arf, 
which sequesters more of Mdm2 and hence reduces the degradation of p53 which therefore 
increases the amount of phosphorylated p53. 
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Atm 




DNA 




Switch 




Damage 




Figure 1: The core regulatory network. The three modules labeled as A, B and C, refer to 
the DNA damage repair, the Atm kinase switch, and the p53 oscillator respectively. These 
are described in detail in the text. Panel l^Spresents a cartoon of the regulation for the 
integrated p53 network. 
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Figure 2: Bifurcation diagram and time series plot for the Atm kinase activity exhibiting 
bistabihty. The upper plot shows the steady state values of Atm-P as a function of the 
DNA damage signal Sd- The thick lines indicate the stable points, whereas the dotted lines 
indicate the unstable points. For Sd > 2.7, the system switches from a low value of Atm-P, 
to a high value of Atm-P (a similar bistable effect is observed for Nbsl-P). The lower plot 
shows the behavior of Atm-P as a function of time. In this simulation, Sn was modelled to 
be an exponentially decreasing function of time, and hence used as a parameter in Eq. 2. 
Initially the damage causes Atm-P to switch to a high level (~ 1), and then subsequently, 
as Sd decreases and crosses a certain threshold, the Atm-P switches back to a low value. 
The plot shows that as damage crosses ^ 0.08, the switch is turned OFF, and one obtains 
a sharp crossover. 
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Figure 3: p53 and Mdm2 dynamics. Panel A shows the steady state value of p53-P as a 
function of Atm-P which shows a supercritical Hopf bifurcation (indicated by point HB) that 
occurs at Atm-P ~ 0.22. At this point the system develops stable oscillations, where the 
amplitudes of the oscillations increase with increasing Atm-P. The limits of the oscillations 
are marked by the thick lines. The dashed line after the point marked HB indicates the 
unstable branch. In panel B, we display the time series of p53-P and Mdm2, which can be 
seen to be out of phase. The overlap between the time series, occurs due to the activation- 
inhibition loop described in the model. For the plotted time series, the value of Atm-P = 1.0. 
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Figure 4: Effects of DNA breaks on p53 and Mdm2 dynamics. Panel A: Time- series plots 

for the concentrations of time courses for Atm-P and p53, for Sd = 68, at t = 0. Panels 

B-D: Time-series plots for p53-P, Mdm2 and Sd as a function of various initial damage 

levels, Soit = 0) = 68, 25, 10. The number of pulses increases with damage, while retaining 

their amplitude and frequency (in our model, the number of pulses is not a linear function 

of the damage. If however, we assume a Michaelis-Menten type of rate law for the damage 

c 

such as assumed in P2|, i-e the rate of change of damage is, oc -[pSS-P]-^-^^^, then for 

small kdam, the rate of decay of damage ( for Sd> 1) is proportional to the concentration 
of p53-P (since Nbsl-P, Atm-P^ 1). From this it follows that the number of p53-P pulses 
required for repair will vary approximately linearly with respect to the amount of initial 
DNA damage.) 



28 



Sensitivity of Model Parameters 
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Figure 5: Sensitivity of model parameters to the system dynamics. A plot of the range of 
values of several parameters in logarithmic scale (base 10), vs the parameter index. The 
parameters are numbered from 1-20, and are listed below. The range (length of vertical line) 
indicates the set of values for which the system exhibits similar dynamics. The square symbol 
in each vertical line indicates the nominal value of the parameter in the model about which 
the range is computed. A short line segment indicates higher sensitivity to changes in the 
dynamics of the model, whereas a long line segment indicates robustness. For example, the 
model is very robust to changes in km (index 10), which is the Michaelis-Menten constant for 
the rate of production of Mdm2. Since the hill coefficient for this rate is 4, the transcription 
rate is fairly independent of the value of km. The model is sensitive to changes in cq (index 
12), the rate of production of p53. The steady state concentration of p53 is relatively low 
in the cell, and this restricts the value of Cq from assuming large values, whereas extremely 
low values of Cq are insufficient to provide enough p53, p53-P to start the oscillations. The 
sensitivity plots were obtained by a combination of bifurcation studies and simulations for 
each of the following parameters in the same order as the "Parameter Index" in the figure. 
a, ko, /ci, /c2, ks, k^, P2, kg, -Pi, kp, P^, cq, pi, kd, kip, kp^, k4a, kd, kp, k^ 
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Figure 6: Bifurcation diagram for the steady state values of Atm-P as a function of Sjj. Curve 
A: A bifurcation plot, for modified Phosphatase activity, which is obtained by changing kQ. 
The upper branch of the steady state values of Atm-P do not fall to zero, at low values of 
the DNA damage signal Sd, which means that as the damage is corrected the Atm switch 
is still in the ON state, and hence p53 oscillations continue even when damage is corrected. 
Curve B: The steady state values of Atm-P are reasonably higher at low values of the DNA 
damage signal Sd, as compared to a similar plot in Fig. 2, panel A. These higher basal 
rates of Atm-P and Nbsl-P would lead to faster DNA repair, for the same amount of small 
DNA damage, when compared to the case considered in Fig. 2, panel A. This plot has 
been obtained by increasing the values of a, (3 which corresponds to higher basal levels of 
phosphorylation (see supplementary materials for parameter values) 
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Figure 7: In the upper plot is shown the steady state values of Atm-P as a function of 
the DNA damage signal So, which shows ultrasensitive behavior. The lower plot shows 
Atm-P as a function of time. Sd, whose initial value is Sd = 5, is modeled here, as an 
exponentially decreasing function in time. Initially, Atm-P jumps almost immediately to a 
high level (~ 1), and then subsequently, as Sd decreases it reduces to a low value, but not 
rapidly enough, as is required by the integrated model. 




Figure 8: Panel A shows the steady state value of p53-P as a function of Atm-P which shows 
two supercritical Hopf bifurcations (indicated by point HB) that occurs at Atm-P ^ 0.6, 2.0. 
Within this region, the system develops stable oscillations. The limits of the oscillations 
are marked by the thick lines. The dashed line between the points marked HB indicate the 
unstable branch. In panel B, we display the time series of p53-P and Mdm2 , for a value of 
Atm-P = 1.0. 
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